// Copyright (C) Kumo inc. and its affiliates.
// Author: Jeff.li lijippy@163.com
// All rights reserved.
// This program is free software: you can redistribute it and/or modify
// it under the terms of the GNU Affero General Public License as published
// by the Free Software Foundation, either version 3 of the License, or
// (at your option) any later version.
//
// This program is distributed in the hope that it will be useful,
// but WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
// GNU Affero General Public License for more details.
//
// You should have received a copy of the GNU Affero General Public License
// along with this program.  If not, see <https://www.gnu.org/licenses/>.
//


#pragma once

#include <random>

#include <pollux/common/base/exceptions.h>

namespace kumo::pollux::functions {
    /// Zeta distribution, the PMF is given by
    ///
    /// \f[ P(X = k) = \frac{1}{k^s\zeta(s)} \f]
    ///
    /// where k is positive integer and \f$ \zeta(s) \f$ is the Riemann zeta
    /// function
    ///
    /// \f[ \zeta(s) = \sum_{i=1}^n \frac{1}{i^s} \f]
    ///
    /// This is mainly used to generate test data with skewed frequencies.
    struct ZetaDistribution {
        ZetaDistribution(double s, int n) : cdf_(n) {
            POLLUX_CHECK(s > 1 && n >= 1);
            double z = 0;
            for (int i = 1; i <= n; ++i) {
                z += pow(i, -s);
                cdf_[i - 1] = z;
            }
            cdf_.pop_back();
            for (double &p: cdf_) {
                p /= z;
            }
        }

        template<typename Generator>
        int operator()(Generator &g) {
            double z = uniform_(g);
            auto it = std::lower_bound(cdf_.begin(), cdf_.end(), z);
            return 1 + (it - cdf_.begin());
        }

    private:
        std::uniform_real_distribution<> uniform_{0, 1};
        std::vector<double> cdf_;
    };
} // namespace kumo::pollux::functions
